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Abstract 

This paper deals with two domain decomposition methods for two dimensional linear Schrodinger 
equation, the Schwarz waveform relaxation method and the domain decomposition in space method. 

After presenting the classical algorithms, we propose a new algorithm for the Schrodinger equation 
with constant potential and a preconditioned algorithm for the general Schrodinger equation. These 
algorithms are studied numerically. The experiments show that the two new algorithms improve the 
convergence rate and reduce the computation time. Besides the traditional Robin transmission condi¬ 
tion, we also propose to use a newly constructed absorbing condition as the transmission condition. 

Keywords. Schrodinger equation, Schwarz waveform relaxation method, domain decomposition in 
space method 


1 Introduction 


The aim of this paper is to apply domain decomposition algorithms to the two dimensional linear 
Schrodinger equation defined on (0, T) X U with a real potential V(t,x,y ) 


( d£u := ( idt + A + V)u = 0, (f, x, y ) £ (0, T) x U, 
\ u(0, x, y ) = u 0 (x, y), (x, y) <G U, 


where U = ( xi,x r ) X ( y b ,y u ) is a bounded spatial domain of M 2 with xi,x ri y b ,y u G M and the initial 
datum uq £ L 2 (U). The equation is complemented with homogeneous Neumann boundary condition on 
bottom and top boundaries and Fourier-Robin boundary conditions in left and right boundaries. They 
read: 


d n u = 0, y = y b , y u , d n u + S b u = 0, x = x h x r , 
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where d n denotes the normal derivative, n being the outwardly unit vector on the boundary <9D, and the 
operator Sb is some transmission operator. 

We consider in this paper two domain decomposition methods. The first one is the Schwarz waveform 
relaxation method without overlap (SWR) fl 2 l fTOj . which is based on the time-space domain decomposi¬ 
tion. The time-space domain (0, T) X Q is decomposed into some subdomains (0, T) X IX,-, j = 1, 2, 

The solution is computed on each subdomain and the time-space boundary values are transmitted via 
transmission conditions. The derivation of efficient transmission conditions is one of the key points of the 
SWR method. For Schrodinger equation, some transmission conditions are proposed in ED El El, such as 
Robin transmission condition, optimal transmission condition etc.. 

The second method we consider here is the domain decomposition in space method (DDS) [7, 14]. First 
of all, the time dependent equation is semi-discretized in time with an implicit scheme on the entire spatial 
domain. This procedure leads to a stationary equation in space. Then standard domain decomposition 
methods (such as the optimized Schwarz method [ 8 ] ( 6 j [133) are applied to this stationary equation. The 
DDS method demands a conforming time discretization. The use of nonconforming discretization in time 
is non standard for the Schrodinger equation and we therefore fulfill this requirement. 

We propose in this article to show the effectiveness of new transmissions conditions (expressed in 
term of absorbing transparent conditions) when we apply the two classical algorithms to the Schrodinger 
equation. The study of the interface problem allows us to introduce some new algorithms which sig¬ 
nificantly reduce both the computational time and the number of iterations. We compare them to the 
classical widely used Robin transmission condition with various intensive numerical tests made on parallel 
computers with up to 1024 subdomains. 

This paper is organized as follows. In Section [2j we present the classical SWR and DDS algorithms. 
We show how the classical DDS algorithm can be interpreted as a combination of some classical SWR 
algorithms. In Section [3j we construct an interface problem and analyse its properties. The discretization 
of the Schrodinger equation is also provided. Based on these properties, we propose new algorithms and 
preconditioned algorithms in the two following sections. In Section | 6 j we provide numerical experiments 
which show the efficiency of our new algorithms. A conclusion is drawn in the last section. 

2 Domain decomposition algorithms 

2.1 Geometric configuration 

The interval ( xi,x r ) is divided into N subintervals ( dj,bj ) without overlap. The points dj and bj denote 
the ends of the subintervals (dj,bj). Thus, the entire domain O is decomposed into N non overlapping 
subdomains Qj = ( dj,bj ) x (yb,y u ), j = 1,2, ...,1V (see Figure [l] for N = 3). We denote the normal 
derivative on subdomain £lj by d nj . 



y 










n 2 



n 2 




SI 2 


C 3 







X 

XI 

II 

II 

= a 2 

b 2 = 

= a 3 b 2 = 

= x r 


Figure 1: Geometric configuration. 

There are obviously other ways to decompose the entire domain. One way is illustrated in Figure 
[2] (left) for N = 4. The intervals ( xi,x r ) and ( yb,Uu ) are simultaneously decomposed into subintervals 
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in both spatial directions. In this configuration, an artificial cross point appears. It is well known that 
the domain decomposition method with cross points is a difficult problem since the problem becomes 
singular at this point. Another possibility is illustrated in Figure [2] (right) for N = 3. The entire 
domain is decomposed into an ellipsis and some rings. This approach has many disadvantages for parallel 
computing. Indeed, we would like to control the number of cells for the meshes of each subdomains. Their 
sizes have to be equivalent to insure a good balance between different process. Thus, we restrict ourselves 
in this paper to the first description (see Figure |TJ) . 
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Figure 2: Two other ways of domain decomposition. 


2.2 Classical SWR algorithm 


The Schwarz methods being iterative, we label the iteration number of the algorithms by k and denotes 
by Uj the solution on subdomain (0, T ) x Qj at iteration k = 1, 2,.... 

The classical SWR algorithm is given by 


( 2 ) 


' 2Sf u k = 0, 

Uj(0,x,y) = u 0 (x,y), 

< BjUj = BjU , 

Bju) = 

. d ^ u j = °> 


{t,x,y) <E (0, T) x Uj, 

(x,y) e ttj, 

x = aj, 
x = bj, 

y = yi,yb, 


with a special treatment for the two extreme subdomains (0, T) x tR and (0, T) x since the boundary 
conditions are imposed on x = ai and x = 6/v 


Biu\ = (<9 n + S b )u\ = 0, x = ai, B N u k N = (<9 n + S b )u k N = 0, x = b N . 


The boundary condition at interface nodes a,j and bj is given in term of operator Bj defined by 

( 3 ) Bj = < 9 n . + Sj, j = 1 , ..., N, 

where Sj is a transmission operator. Various Sj operators can be considered. The classical widely used 
Robin transmission condition is given by 

(4) Sj = -ip, p 6 M + , j = 1,2,..., IV. 

Traditionally, the optimal transmission operator is given in term of transparent boundary conditions. For 
the general linear two dimensional Schrodinger equation, we only have access to approximated version of 
the TBCs given by the recently constructed absorbing boundary condition S™ ade |1] |2j which we used as 
the transmission condition 

(5) Sj = -iy/id t + A Fj + V, j = 1, 2,..., N, 
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where Ti = {61} x (y&,y u ), Ij = {uj,bj} x ( yb,Vu),j = 2, 3, iV - 1 and T^ = {oat} X {yb,y u )- In 
our case, the Laplace-Beltrami operator Ap^. is <9 2 . Numerically, this operator is approximated by Pade 
approximation of order m 


yjidt + A r . + Vu ~ ( V»i'' • V a?d?(idt + A r . + V + 

s =0 s=l 


If the potential V = 0 and the spatial domain is M 2 , the absorbing boundary condition is actually a 
transparent boundary conditions. Then it could be proven that the transmission condition leads to 

an optimal SWR method. 

The classical SWR algorithm is initialized by an initial guess of BjU®\ x=aji b jf j = 1,2,...,A. The 
boundary conditions for any subdomain f lj at iteration k + 1 involve the knowledge of the values of the 
functions on adjacent subdomains 12j_i and at prior iteration k. Thanks to the initial guess, we can 
solve the Schrodinger equation on each subdomain, allowing to build the new boundary conditions for 
the next step, communicating them to other subdomains. This procedure is summarized in (J6| for A = 3 
subdomains at iteration k. 


( 6 ) 


( Biu’H^bA 

B2U2 \x=a 2 

B2U2 I £=62 

\B 3 U 3 \x=a, 3 / 


Solve 



Let us define the flux at iteration k by 


(B 2 u\ u =6l \ 
B\V,2 I £=<22 
B 3 U 2 \x=b 2 
\B2U^\ x= a 3 / 


Comm^ 


/Rr^ +1 U =6l \ 

B 2 U k 2 +1 \ x= a 2 

B2U2 +1 \x=b 2 

\B 3 u k 3 +1 \ x= J 


9 — (#l%U=6i, • • • > Bj Uj | x=aj 5 Bj Uj | x=bj ) ' 1 B]\fUj^\x=a]y') j 


"• T " denoting the transpose of the vector. Thanks to this definition, we give a new interpretation to the 
algorithm which can be written as 


(7) 


9 k+1 = K c g k , 


where 1Z C is a linear operator. The solution to this iteration process is given as the solution to the 
continuous interface problem 

(8) (/ - TZ c )g = 0, g = lim g k . 

k —^00 

where I is identity operator. 


2.3 Classical DDS algorithm 

The other algorithm we consider in this paper is the domain decomposition in space algorithm (DDS). 
The equation ([!]) is first semi-discretized in time on the entire domain (0, T ) x 12. The time interval (0, T) 
is discretized uniformly with Ap intervals of length At. We denotes u n (resp. V n ) an approximation of 
the solution u (resp. V ) at time t n = nAt. The Crank-Nicolson scheme on (0, T) x 12 reads 


.u n 
1 — 


Un —1 

At 


T u n — 1 V n T Vn— 1 u n T u n — 1 

+ 9 9 9 


0, 1 ^ n ^ Ap. 


By introducing new variables v n = (u n + u n -\)/2 with Vo = uq and W n = (V n + V n ~ i)/2, we get a 
stationary equation defined on 12 with unknown v n 


(9) 


r _ 2^ 

b'xVn — 1 ’ 
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where £ x := ^ + A + W n . We recover the original unknown by u n = 2v n — u n -\. Then, the optimized 
Schwarz algorithm is applied to the stationary equation ©• We denote by Rj, j = 1, 2 ,N the restriction 
operator from 0 to ftp At time t n , the classical algorithm reads 


( 10 ) 




°^ >xV n,j ~ At Un ~ 1 > 

j v n,j + Sj v n,j = ^ n j V n,j -1 

A V n,j + SjVnj = ^Aj+l 
d *A,J = °’ 


+ Sjv 

+ SjV 


k -1 

71,7 — 1 ’ 

fc—1 
71,J + 1 ’ 


(x,y) £ %, 
X = CLj, 
x = bj, 

y = yuVb, 


where v k j denotes the unknown at time t n , on subdomain Qj at iteration k and Sj is the semi-discrete 
form of Sj given by Q or (|5]). A special treatment for the two extreme subdomains is needed and the 
boundary conditions read d n] v k ] + S\v k 1 = 0, x = a±, d nN v k N + N = 0, x = b n- 

Since the interval (t n -i,t n ) contains only one time step, the DDS algorithm can be numerically 
interpreted as a sequence of some SWR algorithms. 


Algorithm 1 : DDS algorithm 
The initial datum is uq. 
for 77. = 1 , 2 , do 

Apply the SWR algorithm to 

f Sfu = 0, (t, x,y) £ (t n -i,t n ) x Q, 
\ u(0,x,y) = u n -i(x,y), (x,y)eQ, 


where t n = nAt. 


3 Discrete interface problem 

We know by ([8]) that the classical SWR algorithm reduces to the interface problem 

(/ - lZ c )g = 0, g = lim g k . 

k —^oo 

The aim of this section is discretize this relation and to show that the discrete interface problem can be 
written as 


( 11 ) (I~£h) g = d, 

where the vector g is the discrete form of g, d is a vector 

d = (d J r , d T , dj r , • • • , d^) T £ K NyXNr , 
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and Ch is a block matrix (the notation “MPI j” above the columns of the matrix will be used in section 

0 


( 12 ) 


Ch = 


/MPI 0 

MP] 

1 

MPI 

2 

MPI 

IV-2 

MPI IV-1\ 


x 2 ’ 1 

x 2 ’ 2 






X 1 - 4 











x 3 ’ 1 

x 3 ’ 2 





X 2 ’ 3 

X 2 ’ 4 









X 3 ’ 3 

X 3 - 4 









X^v-i’i 

X N~ 1.2 









X N.l 

V 





* 

CO 

IV—1,4 

/ 


The concrete definitions of d ji, d hr . j = 1,2,...,IV and the bloks in Ch are given in proposition 3.5 


and proposition 3.6 The subsection 3.2 is devoted to the derivation of some important properties of Chi 
which give the ground to make the new algorithm given in section [4] 


3.1 Preliminaries related to discretization 

Without loss of generality, we present the discretization of Q on (0, T) xH with the following boundary 

conditions mm 

d n u = 0, y = ybiVui 
d n u + Su = l(t, y ), x = xi, 
d n ii + Su = r(t , y)i x = x r , 

where l(t,y ) and r(t,y ) are two functions. The discrete version of ([2]) follows immediately. Concerning 
the semi-discretization with respect to time, we follow the procedure given in section [2~3| for the classical 
DDS algorithm. The scheme therefore reads 


(13) 


2 i 2 i 


The semi-discrete transmission condition is given by 

^nVn T Si) yi — l n i X — Xh dnV n T Sv n — V n i X — X ri 

where l n = l{nAt,y), r n = r(nAt,y ) and S is the semi-discrete form of S. If we consider the Robin 
transmission condition Q, we have 

(14) Robin : Sv n = —ip ■ v n . 

The approximation of the transmission condition ([5j is given by 


(15) 


cm 

*^pade ' 


sv n = -% y <v n +* y «r«~ 1/2 , 

s=0 s=l 


where a™ = e l9 / 2 /(m cos 2 (^^^)), d™ = e^tan 2 ^ 2 ^^ ), s = 0, 2,..., m, 9 = |. The auxiliary functions 
p’s 1//2 , s = 1,2,..., m are defined as the solutions of the set of equations 

(^ + Ar + W n + C)^r 1/2 -v n = ^r\ 


•At 

„n r) n —1/2 ,„n— 1 ,„0 n 
Vs = 2 Vs Vs , Vs = 
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The spatial approximation is realized by the standard <Qh finite element method. The uniform mesh 
size of a discrete element is (Ax, Ay). We denote by N x (resp. N y ) the number of nodes in x (resp. y) 
direction on each subdomain. Let us denote by v n (resp. u n ) the nodal interpolation vector of v n (resp. 
u n ), l n (resp. r n ) the nodal interpolation vector of l n (resp. r n ), M the mass matrix, S the stiffness 
matrix and Mw n the generalized mass matrix with respect to W n v(f>dx. Let M r the boundary mass 
matrix, § r the boundary stiffness matrix and ML the generalized boundary mass matrix with respect 
to / r W n V(pdT. We denote by Qi (resp. Q r ) the restriction operators (matrix) from 17 to {x;} X ( ybiVu ) 
(resp. {x r } x (j/6, y u )) and Q T = ( Qj , Qj). The matrix formulation for the transmission condition Robin 
is therefore given by 


(16) 


Robin : 


A n + ip ■ 


2 i 

v n = —Mu n _i - . 


I r Q T 


where A n = — S + Mw n . The size of this linear system is N x X N y . If we consider the transmission 

condition <S'™ de , we have 


(17) 


with 


/A n + i{J2T=o a 

C 

c 

c 


QM l Q 


rnT 


QM l Q 


rnT 


QM r Q T J 


A 


/ U n _l\ 

,jn—x 

<P1 

. _n—1 

F>2 


Wm 1 / 


/ v n \ 

n— 1/2 

v i 

n— 1/2 
¥>2 


/ W 1/2 / 


( M r Q T 

0 

0 


B s = —iaf^ dff'MF Q T , 1 ^ s ^ m, 

C = -QM r , 

= Q(|^M r - §r + M w n + d™M r )Q T , 1 < s < m. 

It is a linear system with unknown (v n , ip™ 1/<2 ,..., 1 ^ 2 ) where 1 ^ 2 is the nodal interpolation of 

ifs on the boundary. 

Remark 3.1. The S'™ c | e transmission condition involves a larger linear system to solve than the one of 
the Robin transmission condition. The cost of the algorithm with the 5™^ transmission condition is 
therefore more expensive. 

The equations (16) and © are given for a fixed discrete time t n . They however can be written 
globally in time by some straightforward calculations. 


Proposition 3.2. For the Robin transmission condition, the global form in time of the equation (16) is 
(18) (A - B)v = F - M r Q g, 
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where B = -ip • M 1 = -ip • diag^jM 1 }, Q( = diag Nt {Q{ }, = diag^jQ,! } and 


( Ai 

-S M A 2 


A = 


4 i ■ 
At 1 




Aa 


( v i \ 
Vo 


V = 



V 


-f t m a Nt 


\vn t J 



/ Muo \ 


(QJ 


Qi 

\ 

At 

— Muo 

, Q t = 

Qi 


Qi 





V 

Qi 


Qi) 


Proposition 3.3. If we consider the transmission condition S™: vie , then the equation ( |17[ ) can be written 
globally in time as (18) 

with B given by 


(19) 


B = - 




(A - B)v = 

F-M r Q T g, 


l C 0 

M r + Y 1 ’ 1 



\ 


Y 2 ’ 1 

c a M r + Y 2 ’ 2 




Y 3 ’ 1 

Y 3 ’ 1 

c a M r + Y 3,3 


V 

Y-Wh 

Y N T,2 

Y^T;3 

• c a M r + Y Nt ’ Nt j 


where c a = i{Y!T=o a T)- 

Proof. According to ( | 17 [ ), for s = 1 , 2 , n = 1 , 2 , Nt, we have 

-QM r v n + D>y - ;1/2 = ^QjM r Q T ip ™ -1 

where B” = Q(^M r — § r + + d™M r )< 5 T . We therefore obtain 

2 i 


^n- 1/2 = m -l QjM r Vn + 


and 


<P. 


r 1 = vr 3/2 - <^r 2 , *>2 = o, n > 2 . 


^r 1/2 = E L " ,p v^ 

p=i 


By induction, ys”' is given by 

( 20 ) 

where lff’ p are matrix. Replacing tpf ^ 2 by v p in the first row of the equation © gives 
(21) (A n + i(£ a™) • M r )v n + £ Y n ’ p v p = ^Mu n _! - M r Q T (”) , 

s =0 p= 1 ^ ' 

where Y n,p := —z a™d’J l M r Q r L"’ p Then, according to ( [2l| ) , the matrix B is given by (19). 

Remark 3.4. Following the same procedure, the equation ([2]) is discretized on each subdomain (0, T ) X f lj. 
Accordingly, we define the matrices A j, B ; , M r T Lj’j associated with the hnite element method, the 
restriction matrix 0,3-1- Q j,r and the solution vector v^ n . The subscript "j" emphasizes the definition of 
the matrices associated to the subdomains (0, T) X f lj. 
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3.2 Properties of Ch 

Let us define the fluxes 


l*(t, y) = d nj v k + SjVj , rj(t, y) = d nj v k + V/> (*, 2/) € (0, T) x 

for j = 1, 2,iV with the two special cases l\ = = 0 corresponding to the left and right boundaries 

and 

ln,j = lj(tn,y), = rj(t n , y). 

We denote by ft n (resp. ) the nodal interpolation vector of ft (resp. r k n ) 


j,n \ t 3i n/ * - j,n v Jr j,ni 

/C+l _ 1 k | 0/0 . O 


( 22 ) 


■Ti,n = -l£n + 2Q;,z • j = 2,3, N, 


lj + l,n — r j,n + 2Qj,r ' Sj v j. n , j — P 2, IV 1, 


where we have for the Robin transmission condition 


Robin : Sjv| n = -ip • v£ n , 


and for the S™ de condition 


s™,d. : S;v‘„ = -i(£ + j <C<V 


n—1/2 

j,s 


s=0 

Then, we define the discrete interface vector g by 
(23) 
where 


S=1 


g = lim g fc , 
k —^oo 


gN 


g 

k, T 
r l,l > 

- O fc - T 


= ( g r 


fc _ ofc.T _ 

g J W,1 ’ 


g^r = (1^,1, i; 

ife,T Jc,T 


k, T fc.T fc.T\T 

g 2 ,-,gjv ) > 

■,r{J T ) T 6CW, 

ifc,T \T c (r^NyXNx 

> hv,iv T ! ^ > 


fc,T 

1,2 > ' 

fc,T 
TV,2’ ' 


fc,T 

r j,Nx 


) T e C 2N * xNt , j = 2,3, ...N — 1. 


Proposition 3.5. If we consider the Robin transmission condition, the discrete form of the interface 
problem ([8]) is given by & 

Proof. According to (|22|, and the definitions of g fc , it is easy to verify that 


(24) 

X jA = -I- 2 ip • Qj,i(Aj - B/) 'M' Q;,. 

A J ’ 2 = — 2ip • Qj,i(Aj - By) 1 M 1 

X j ’ 3 = -2ip • Qj, r (Aj - Bj) _1 M r ^QT,, 

X jA = -/ - 2ip ■ - B,) 'M'df ,.. 


and 


(25) 

dyi = 2ip • Qj-i, r (Aj_i - Bj—1) —1 Fj_i, j = 2,3, ..., N, 
d j t r — 2ip • Qj+i i /(Aj_)_i Bj^i) F,/ ■ i ■ j — 1, 2, ..., N 1 
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Proposition 3.6. If we consider the transmission condition 5™^, then ( |11[ ) is the discrete form of the 
interface problem (pi). 


Proof. By using (22), we have 


Sp,.n = -C„v‘„ + i £ <C<C £ L”fv‘„ = +, £ (£ )r 

p =1 P=1 S=1 


k 

j,p ' 


S=1 


We could easily verify that 


(26) 


X*’ 1 = -I + 2Q^Bf(A j - B,) 

X-t’ 2 = 2Q jii Bf (A j - Bj) _1 M rj Qj r , 

A''' 3 = 2Q,,.ByiA, - lij) i M r 'Q / : /: 

W’ 4 = -/ + 2Q ;; I ij ( A, - B/) 'M r Q ;> . 


and 

(27) 


dj,i = - B i-i) 'F.-i, j = 2, 3,A, 

d j,r = Qj+i,/B| +1 (Aj + i — Bj + i) _1 Fj + i, j = 1, 2, N — 1, 

where the matrix Bj is defined by 

( 


Bj> = 


-Cal + i^aTd? m 

S = 1 
m 


-C, 


. a i+iJ2 a T^ L ?; s 2 


S=1 


S=1 






A^t,1 




Nt,2 

‘j,s 


-cJ + i^aTd™ I* 


' m.jr Nt,Nt 


s=l 


S=1 


S=1 


/ 


Let us now study the structure of the subblock of Ch for time independent potential V = V(x,y). 
More specifically, we focus on 


(28) 


A 1 ’ 4 = {^’ 4 }l<n, S <AT T , 

X jA = {x{d s }l^ n ,s^N T ,X j ’ 2 = {x{; 2 s }i^ n ,s^N T , 

x j ’ 3 = (4’fj i^ Nt ,X^ = K’fji ^N T ,j = 2,3,..., N — 1, 

vJV,l _ ; JV,li 

x*- \^n,s jl • 


where xfl] s , xf^ s . xiif s , xii s £ C NyXNy are submatrices. 


Each subblock Ay 


1 , 24,4 


are made of submatrices which are set on very specific positions. This structure 
is presented in Figure [3] for 3 time steps and 6 nodes on the interface between two subdomains. We see 


that each sub-diagonal block is identical. We present this property mathematically in proposition 3.7 


with the two transmission condition Robin or S™ de . The demonstration is similar to the one obtained for 
one dimensional Schrodinger equation [5], The formal difference between dimension one and dimension 
two is that the flux are scalar in one dimension and vectors in two dimensions. 
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6x6 



6x6 

6x6 


6x6 

6x6 

6x6 


Figure 3: Block structure, Nt = 3, N y = 6. 


Proposition 3.7. For the transmission Robin (resp. S™ ade ), assuming that the linear system (16) (resp. 

is not singular, if the potential is time independent V = V(x,y), then the matrices A' 1,4 , X^ 1 ,X^ 2 , 
X ■ ?,3 ; X J ’ 4 , j = 2,3, ...,N — 1 and X N}1 are block lower triangular matrices and they satisfy 


„1,4 


1,4 


7,1 7,1 1,2 7,2 

rytj J - /'V**' 7 7y*i/ 7 - 7 

J 'n,s J 'ri-l,s-l> 


r j,3 _ J,3 j,4 _ j,* 7=23 AT — 1 

b n,s ,b n— l,s—1’ x n,s h n— l,s— 1>J .l, 

7V,1 _ AT,1 
■^n,s ~ ®ra—l,s—1 > 

/or 2 ^ s ^ n ^ IVt- 

If the potential is constant, the sub-blocks of Ch are simpler and follow the following proposition. 

Proposition 3.8. Assuming that the matrix A j — B j is not singular with the Robin or the S™ ade trans¬ 
mission condition, and that the mesh is uniform and the size of subdomains Qj are equal, if the potential 
V is constant, then the subblocks of Ch satisfy 


n— 1,3—1’ 
j,4 _ jA 


(29) 


j^-2,1 _ jj^-3,1 _ . _ X^^ X^^ _ X^^ _ . _ j^N—1,2 

j^-2,3 _ x 3,3 _ ... _ j^N-1,3 j£-1A _ x _ ... _ x^ —]-A 


Proof. Thanks to the hypothesis of the proposition, the geometry of each subdomain is identical. Thus, 
the various matrices coming from the assembly of the finite element methods are the same. Therefore, 
we have 

Ai = A 2 = ■ ■ ■ = Ajv, M ri = M r2 = ■ ■ ■ = M FiV , 
and the restrictions matrices satisfy 

Ql,l = Q2,l = • • • = Qn,U Ql,r = Q2,r = • • • = Qn,t- 
Since V is constant, then M]j/ n = M 2 j u/„ = ••• = Thus by definitions, we have 

Ai = A 2 = • • • = Atv, Bi = B 2 = • • • = Bat, Qi = Q 2 = • • • = Qyv, 


and Bf = B| = • • • = B^ for the transmission condition S'' 


(24) and (26). 


pade* 


The conclusion therefore follows from 
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4 New algorithms for the Schrodinger equation with a constant poten¬ 
tial 


Thanks to the analysis yielded in previous section, we can build explicitly the interface problem 
The main idea of our new algorithm is therefore to explicitly construct the matrix Ch and the vector d. 
Following propositions (3.7) and (3.8), it is sufficient to compute four subblocks to explicitly build the 
matrix Ch- Without loss of generality, we construct the blocks X 2,1 , X 2,2 , X 2,3 and X 2,4 . Furthermore, 
only the first N y columns of each block are necessarily computed. 


4.1 New SWR algorithm 

We propose here a new SWR algorithm for the Schrodinger equation with a constant potential V. We 
always suppose that the size of each subdomain is identical and the mesh is uniform. 

Algorithm 2: New SWR algorithm, constant potential V 
1 : Build Ch and d in relation (/ — Ch )g = d explicitly. 

2 : Solve the linear system (/ — Ch) g = d by an iterative method (fixed point method or 
Krylov methods). 

3: Solve the Schrodinger equation on each subdomain on (0, T) X Qj using the flux 
obtained from step 2. 


Mathematically, this new SWR algorithm is identical to the classical one. The main novelty here is that we 
construct explicitly the matrix Ch and the vector d (corresponding to the first step), while in the classical 
algorithm, Ch remains an abstract operator. It is not usual to build explicitly such huge operator, but as 
we will see, its computation is not costly. 

We show below the construction of Ch and d for the Robin transmission condition. This is is based 
on the formulas (24) and (25). For the transmission condition S , ™ de , the idea is similar, but involves (26) 
and (27). According to the proposition 3.5 the column s of X 2,1 and A' 2,3 are 


X 2,1 e s = —e s -2 ip- Q 2> ;M r2 Q^M r2 Q^e s , 
X 2 ’ 3 e s = —2 ip ■ Q 2 ,r(A 2 - B 2 ) _1 M r2 Q^e s , 


where the vector e s = (0, 0,..., 1, ...0) E (^N T xN y ^ a jj q s elements being zero to the exception of the s-th, 
which is one. The element M lz Q^e s being a vector, it is necessary to compute one time the application 
of (A 2 - Ba)" 1 to a vector. Similarly, we have 


X 2,2 e s = -2 ip ■ Q 2 ,*(A 2 - B 2 ) _1 M r2 Q^ r e s , 

X 2,4 e s = —e s - 2 ip ■ Q 2 , r (A 2 - B 2 ) -1 M r2 Q^ r e s , 

In conclusion, to know the first N y columns of X 2,1 , X 2,2 , X 2,3 and X 2,4 , it is sufficient to compute 2 N y 
times the application of (A 2 — B 2 ) _1 to a vector. In other words, this amounts to solve the Schrodinger 
equation on a single subdomain 2 N y times to build the matrix Ch- Without loss of generality, since the 
geometry of each sub domain is identical and the matrix Ch is independent of the initial solution, we focus 
here on the domain (0,T) X 0 2 . The resolutions being all independent, we can solve them on different 
processors using MPI paradigm. We fix one MPI process per domain. To construct the matrix Ch , we 
use the N MPI processes to solve the equation on a single subdomain ((0, T ) x 0 2 ) 2 N y times. Each MPI 
process therefore solves the Schrodinger equation on a single subdomain maximum 


2 N„ 


X mp i ] 


+ 1 times, 
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where [x] denotes the integer part of x. This construction is therefore super-scalable in theory. Indeed, if 
N is doubled, then the size of the subdomains is divided by two and A mp i is also approximately halved. 

Concerning the computational phase, the transposed matrix of Ch is stored in a distributed manner 
using the PETSc library |4j. As shown by (12), the first block column of Ch lies in MPI process 0. The 


second and third blocks columns are in MPI process 1, and so on for other processes. The size of each 
block is (Nt x N y ) x (Nt x N y ). Each block contains (Nt + 1) x Nt/2 x N y nonzero elements. 

The construction of the vector d is similar. According to (25), one needs to apply (Ay — By ) -1 to 
vector Fj for j = 1, 2, In other words, we solve the Schrodinger equation on each subdomain one 

time. Again, the vector is stored in a distributed manner using the PETSc library. 


Ul A }MPI 0, N t x N y elements 


d = 2 


d i,i 

d j,r 

\&N,l) 


MPI j. 2 Nt X N y elements 


}MPI N , Nt x N y elements 


4.2 New DDS algorithm 

Since we have interpreted the DDS method as sequence of some SWR methods in Algorithm [l] we can 
apply directly the ideas developped in the previous sections to the DDS method. Let us denote the 
interface problem of 

( id t u + Au + Vu = 0, (t,x,y) € (t n -i,t n ) x Q, 

\ u(0,x,y) = u n -i(x,y), (x,y)€tl, 

by 


(30) 


(d Ch t n)Sn — d n , 


where Ch,n is interface matrix and d n is vector. In the classical algorithm, C} lri remains an abstract 
operator. We propose to build it explicitly in the new DDS algorithm. Actually, thanks to the following 
proposition, the computation of the complete operator Ch,n only require to compute Ch, i- 

Proposition 4.1. For the transmission Robin and 5'™ de , the interface matrix satisfies 

Ch,l Ch,2 ■■■ Ch'Nqp' 


Proof. According to (24) and (26), the interface matrix is independent of the initial datum, thus the 
conclusion follows directly. 


Algorithm 3: New DDS algorithm, constant potential V 
1 : Build Ch,i explicitly. 

2 : The initial datum is uq. 
for n = 1,2, ...,Nt do 

2 .1: Build d n on time t n , 

2 .2: Solve the linear system (I — Ch, n )Sh,n = d n by an iterative method, where 
the initial vector is chosen as ghn-l- 

2.3: Solve the Schrodinger equation on each subdomain (t n -i,t n ) X using the 
flux from step 2.2 to compute u n . 
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Mathematically, this new DDS algorithm is identical to the classical one. Compared with the new SWR 
algorithm, the construction of the interface matrix is less costly since the size of Ch.i is smaller than that 
of C h . 


5 Preconditioned algorithms for general linear potential 


In the case of a non constant potential, the proposition |3.8| does not hold. Thus we could not construct 
easily the matrix Ch- The aim of this section is to present preconditioned algorithms for the Schrodinger 
equation (jlj with a non constant potential V = V(t,x,y). Adding a preconditioner to the equation (11) 
leads to a preconditioned SWR algorithm 


(31) 


p- 1 (i-c h ) g = p~ 1 d, 


We propose here to use 


P = I-C 0 , 

where Co denotes the interface matrix in ( [TT| ) defined for the free Schrodinger equation when V = (1. 
As mentioned in the previous section, it could be easily constructed numerically and it is stored in a 
distributed manner. For any vector y, the vector x := P~ 1 y is computed by solving the linear system 


Px = y. 


with Krylov methods (for example GMRES or BiCGStab). However, the size of P increases linearly with 
the number of subdomains N. Also the number of involved operations for multiplying Co and a vector 
(which is the basic operation of Krylov method) is not negligible compared to the computation of the 
solutions to the equations on each subdomain if N is large. Thus, the application of the preconditioner 
increases. 

We could derive straightforwardly a preconditioned DDS algorithm from the point of view of Algorithm 
[l] The preconditioner for all time steps n = 1,2,..., Nt is always chosen as 


P = I-C 0) i, 


where £q,i is the interface matrix of (30) defined for V 


0 and n = 1. 


Remark 5.1. The idea of these new algorithm are not limited to the Schrodinger equation. It could also 
be applied to some other PDEs, like heat equation etc. The only limitation is that the mesh and the 
decomposition should be uniform. 


6 Numerical results 

The complete domain D = (—16,16) x (—8,8) is decomposed into N equal subdomains. The size of a 
single cell is Ax X Ay. We consider two different meshes 

Ax = 1/128, Ay = 1/8; 

Ax = 1/2048, Ay = 1/128. 

With the first mesh, it is possible to solve the Schrodinger equation (JI]) on the entire domain on a single 
node of a cluster composed of 92 nodes (16 cores/node, Intel Sandy Bridge E5-2670, 32GB memory/node). 
Thus we could observe if the parallel algorithms allow to reduce the total computation time of the 
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sequential algorithm. We are interested in the strong scalability up to 1024 subdomains. The initial 
datum in this section is 

(32) u 0 (x,y) = e~ x2 - y2 -°- 5lx . (see Figure [4]) 


This section is composed of two subsections. The first one is devoted to the free Schrodinger equation 
(V = 0). In the second, we consider the Schrodinger equation with a linear potential V = x 2 + y 1 . 

Remark 6.1. For the DDS algorithm, we consider mostly the convergence properties (for example the 
number of iterations) of the first time step. In other words, we only study the evolution between t° —>• t 1 
to compute v l with the optimized Schwarz method. 

Remark 6.2. The theoretical optimal parameter p in the Robin transmission condition is not at hand for 
us. We then seek for the best parameter numerically for each case. 


6.1 The free Schrodinger equation 


In this part, we first compare the SWR method and the DDS algorithm. The DDS algorithm being 
more efficient, we next compare the classical DDS and the new DDS algorithms. Finally, the influence of 


parameters is studied in Section 6.1.3 


6.1.1 SWR vs. DDS 

We compare the SWR and the DDS methods in the framework of the classical algorithm. The time step 
is set to At = 0.01. The mesh size is Ax = 1/128, Ay = 1/8. Both methods are applied on time-space 
domain (0, T) X D with various final time T. We denote by the number of iterations required to 

obtain convergence of the SWR method, N® ei the number of iterations of the first time step of the DDS 
method and T s (resp. T D ) the total computation time of the SWR method (resp. DDS method). Table 
[T] and Table [2] present the numbers of iterations and the computation times of both methods for N = 2 
and N = 32, where the transmission conditions are S™ ade ,m = 5 and Robin respectively. The initial 
vector is the zero vector and the GMRES method is used on the interface problem. We could see that 
N^ ei increases with the final time T and N- der > N® er . Thus T s > T D . 

In the next subsections, without special statement, we only consider the DDS method since it does 
requires less computation time. 
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Table 1: Number of iterations and computation time of the SWR method and the DDS method for 
N = 2, 32 with the transmission condition S^,^. 


T 


N -- 

= 2 



N = 

32 



NLr 

Ni? er 

rj-tS 

rj~\D 

nL 


rj-iS 

rj-iD 

0.05 

17 

9 

17.6 

12.3 

17 

10 

1.5 

1.9 

0.1 

25 

9 

44.5 

21.5 

25 

10 

3.5 

1.8 

0.15 

30 

9 

78.4 

30.9 

31 

10 

6.4 

2.6 

0.2 

44 

9 

147.9 

40.0 

45 

10 

11.8 

3.4 

0.25 

51 

9 

215.0 

49.3 

52 

10 

16.9 

4.1 

0.3 

55 

9 

271.0 

58.6 

55 

10 

21.2 

4.8 

0.35 

58 

9 

332.1 

67.9 

59 

10 

26.3 

5.6 

0.4 

61 

9 

402.9 

77.2 

62 

10 

32.0 

6.3 

0.45 

64 

9 

474.5 

86.4 

65 

10 

37.6 

7.1 

0.5 

68 

9 

557.6 

96.1 

73 

10 

46.3 

7.8 


Table 2: Number of iterations and computation time of the SWR method and the DDS method for 
N = 2,32 with Robin transmission condition. 


T 


N -- 

= 2 



N = 

32 



NLr 

Ng ei 

rj-iS 

rpD 

nL 


rj-iS 

rj-iD 

0.05 

19 

11 

13.0 

15.9 

20 

11 

1.1 

0.6 

0.1 

28 

11 

32.6 

16.8 

29 

11 

2.7 

0.6 

0.15 

44 

11 

73.3 

17.9 

47 

11 

6.1 

0.7 

0.2 

57 

11 

123.1 

18.9 

57 

11 

9.6 

0.8 

0.25 

76 

11 

204.6 

19.9 

80 

11 

16.6 

0.9 

0.3 

85 

11 

271.4 

21.0 

89 

11 

22.1 

1.0 

0.35 

90 

11 

333.6 

21.9 

92 

11 

26.7 

1.0 

0.4 

100 

11 

425.3 

23.0 

102 

11 

33.5 

1.2 

0.45 

109 

11 

519.4 

24.0 

208 

11 

40.0 

1.3 

0.5 

114 

11 

601.2 

25.0 

116 

11 

47.3 

1.3 
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6.1.2 Comparison of classical and new algorithms 

In this part, we are interested in the performance (number of iterations and computation time) of the 
classical and the new algorithms with the two transmission conditions. We observe the strong scalability 
of the two algorithms. Both algorithms and transmission conditions are compared in the framework of 
the DDS method. The final time is T = 0.5 and the time step is fixed as At = 0.01. The GMRES method 
is used on the interface problem. The initial vector is the zero vector. Since there is no theoretical 
result for us on the choice of the optimal parameter p in the Robin transmission condition, we make 
tests with different p to find the numerically optimal one. However, we will see in the next subsection 
that the number of iterations and the computing time are not sensitive to p using the GMRES method 
on the interface problem. Thus, it is difficult to choose the optimal one. We take p = 15 for the mesh 
Ax = 1/128, Ay = 1/8 and p = 10 for the mesh Ax = 1/2048, Ay = 1/64. For the transmission condition 
S'™de, we take m = 5, the numerical optimal according to some tests. 

We first consider the mesh Ax = 1/128, Ay = 1/8. Figure [5] presents the convergence history of the 
first time step for N = 2 and N = 32. We also show the number of iterations of the first time step and the 
total computation time in Table[3]for N = 2,4,8,16, 32. The boundary conditions involve the operator Sj, 
which can be the usual Robin transmission operator or S'™ de . The reference computation times for solving 
the Schrodinger equation ([!]) on the entire domain are therefore not the same. We denote respectively 
“Robin ref” and “S'p ade ref.” the solution of the Schrodinger equation computed on the whole domain on 
a single processor for the Robin transmission condition and for the transmission condition 5'/ ade . 

We see that the new algorithm allows us to reduce the computation time compared with the classical 
algorithm. Moreover, the number of iterations is almost independent of the number of subdomains and 
the algorithm is scalable. Finally, the algorithm converges faster with the transmission condition 5'™ de , 
but takes more computational time since. Indeed, the application of the transmission condition S™ ade is 
more expensive than the Robin transmission condition. 




Figure 5: Convergence history of the first time step for N = 2 (left) ,N = 32 (right). 

Remark 6.3. From the results shown in the following subsection (see Table [5] and Table [6] together), we 
can see clearly the different convergence rates with the two kinds of transmission conditions. 

Next we make tests with the mesh Ax = 1/2048 Ay = 1/64. The entire domain is divided into 
N = 128,256,512,1024 subdomains. We present in Figure [6] the convergence history for N = 256,1024 
and in Table [4] the numbers of iterations and the total computation time. We can see that the classical 
and the new algorithms are not very scalable since the number of iterations increases with the number of 
subdomains. However, the new algorithm takes less computation time. Besides, the number of iterations 
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Table 3: Number of iterations of the first time step and total computation time in seconds, T = 0.5, 
At = 0.01, Ax = 1/128, Ay = 1/8._ 


N 


2 

4 

8 

16 

32 

Number of iterations 

Robin* 

11 

11 

11 

11 

11 



‘-'pade 

9 

9 

9 

9 

10 


Robin ref. 

16.0 


Robin els. 

63.1 

32.6 

17.5 

11.0 

5.4 

Computation time 

Robin new 

30.0 

9.7 

4.4 

2.5 

1.3 

Spa,de ref. 
Side els. 
Spade new 

22.1 


96.1 

49.8 

26.7 

15.0 

7.8 


38.2 

14.7 

6.6 

3.4 

1.8 


required for the transmission condition is less than the one for the Robin transmission condition. 

The computation times are however similar. 




Figure 6: Convergence history of the first time step for N = 256 (left), N = 1024 (right). 


Table 4: Number of iterations of the first time step and total computation time in seconds, T = 0.5, 
At = 0.01, Ax = 1/2048, Ay = 1/64._ 


N 


128 

256 

512 

1024 


Robin* 

14 

16 

22 

36 

Number of iterations 




o5 

‘-'pade 

12 

14 

19 

29 


Robin els. 

250.2 

143.8 

92.5 

101.4 

Computation time 

Robin new 

59.1 

38.1 

36.2 

52.3 


‘-’pade els. 

- 

187.5 

162.6 

127.5 


*^pade H6W 

- 

42.7 

36.2 

45.0 


the memory is insufficient; *: p = 10. 


6.1.3 Influence of parameters 

We study in this part the influence of parameters: m in the transmission condition S/” de and p in the 
transmission condition Robin. The time step is fixed as At = 0.01. The mesh is Ax = 1/128, Ay = 1/8. 
Three different methods are used to solve the the interface problem: fixed point method, GMRES method 
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and BiCGStab method. Two initial vectors are considered: the zero vector and the random vector. In 
our tests, the algorithm initialized with the zero vector converges faster than the one with the random 
vector. However, in this subsection, our goal is to study the influence of parameters. As explained in 
[9], initialization with a zero vector to compute a smooth solution makes that the error contains only low 
frequencies and could therefore draw the wrong conclusions. Thus we consider both of the two initial 
vectors, but we have no wish to compare them. 

The parameter m in the transmission condition S'™ de . The numbers of iterations of the first 
time step with various m are presented in Table [5] We could see that the number of iterations is not 
sensitive to the parameter m if the GMRES method or the BiCGStab method is used on the interface 
problem. However, if the fixed point method is used, increasing the order of S'p" de does not ensure better 
convergence property. The transmission condition S™' de is based on formal Pade approximation of square 
root operator, this approximation may deteriorate for large m. 


Table 5: Number of iterations for different m, N = 32. 



Fixed point 

GMRES 

BiCGStab 

m 

Zero 

Random 

Zero 

Random 

Zero 

Random 

3 

34 

124 

11 

30 

6 

17 

4 

26 

96 

10 

28 

6 

16 

5 

21 

79 

10 

27 

5 

15 

6 

18 

68 

9 

26 

5 

15 

7 

17 

61 

9 

25 

5 

14 

8 

18 

56 

9 

25 

5 

14 

9 

19 

52 

10 

25 

5 

14 

10 

21 

50 

10 

25 

5 

14 

15 

32 

46 

11 

26 

6 

15 

20 

43 

50 

12 

28 

6 
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The parameter p in the transmission condition Robin. We present in Table [6] the numbers of 
iterations with various p for N = 32. From the table, we can see that the algorithm is not sensitive to 
p if the GMRES method or the BiCGStab method is used on the interface problem, while it exists an 
optimal p (marked with an underline) if the fixed point method is used. Besides, the Krylov methods 
could accelerate a lot the convergence. 

In conclusion, the difference between the two transmission conditions is clear if the fixed point method 
is applied to the interface problem and the initial vector is the random vector (see Tables[5]and[6]together). 
The number of iterations for the transmission condition 5'™ dc is less than that for the transmission 
condition Robin. But the difference is smaller using the Krylov methods and the zero vector as the initial 
vector. From the point of view of computation time, the zero vector and the GMRES method are a good 
choice. According to the tests in the previous subsection, the computation time for the two transmission 
conditions are similar. 

6.2 Case of non-zero potential 

The aim of this section is to compare the classical algorithm and the preconditioned algorithm in the 
framework of the DDS method with the fixed point method used on the interface problem. We first 
consider the potential V = x 2 +y 2 . Let us denote by N nopc (resp. N pc ) the number of iterations required 
to obtain convergence of the classical algorithm (resp. the preconditioned algorithm) and T nopc (resp. 
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Table 6: Number of iterations for different p, N = 32. 



Fixed point 

GMRES 

BiCGStab 

p 

Zero 

Random 

Zero 

Random 

Zero 

Random 

5 

57 

580 

11 

35 

6 

21 

10 

34 

315 

11 

32 

6 

19 

15 

32 

239 

11 

31 

6 

18 

20 

35 

209 

11 

31 

6 

19 

25 

40 

200 

11 

32 

6 

19 

30 

46 

199 

11 

32 

6 

19 

35 

52 

204 

12 

33 

6 

19 

40 

59 

209 

12 

33 

6 

20 

45 

65 

222 

12 

34 

6 

20 

15 

26 

32 

198 






T pc ) the computation time of the classical algorithm (resp. the preconditioned algorithm). 

First, we present in Figure[7]the spectral properties of (/ —£) and P^ 1 (I — C) for N = 2, 8 subdomains. 
We see that the eigenvalues of P _1 (I — C) are closer to 1+0 i than those of (/ — £). The convergence 
history is therefore better. The Table [7] shows the number of iterations of the first time step and the 
total computation time to realize a complete simulation. The mesh is Ax = 1/128, Ay =1/8. As 
mentioned before, it is possible to solve the Schrodinger equation on the entire domain with this mesh. 
The computation time is denoted by T rcf . We can see that all algorithms are robust and the number of 
iterations is independent of the number of subdomains. The classical algorithm and the preconditioned 
algorithm are both scalable and the preconditioner allows to reduce the number of iterations and the total 
computation time. In addition, the times of computation T nopc , T pc are less than the reference time T rcf 
if N is large. 




Figure 7: Eigenvalues of the classical (Cls. ) and the preconditioned (Pd.) algorithm, first time step, 
At = 0.01, Ax = 1/128, Ay = 1/8, N = 2 (left), N = 8 (right). 

Next, we consider the mesh Ax = 1/2048, Ay = 1/64. The size of (I — £) is too large to compute all its 
spectral values. The convergence history and the computation time are presented in Figure [8] and Table 
[8] The parameters m used are also presented in Table [8j We can see that the preconditioned algorithm 
is more robust and requires much less number of iterations. However the computation time of the two 
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Table 7: Number of iterations of the first time step and total computation time in seconds of the classical 
algorithm and the preconditioned algorithm, T = 0.5, At = 0.01, Ax = 1/128, Ay = 1/8. 


N 

2 

4 

8 

16 

32 

N nope i TYl = 7 

17 

17 

17 

17 

17 

Npc, m = 5 

5 

5 

5 

5 

5 

rj-i ref 

16.1 

4nopc 

142.7 

75.3 

40.1 

23.9 

12.1 

T P c 

91.3 

43.3 

22.8 

13.1 

7.3 


algorithms are not quite scalable since for the classical algorithm, the number of iterations increases 
with the number of subdomains. Concerning the preconditioned algorithm, the size of preconditioner 
is (2 N — 2) X Nt x N y . This increases with the number of subdomains N. Thus, the application of 
preconditioner takes more computation time with bigger N. 




Figure 8: Convergence history of the first time step of the classical algorithm and the preconditioned 
algorithm, At = 0.01, Ax = 1/2048, Ay = 1/64, N = 256 (left), N = 512 (right). 


Table 8: Computation time in seconds of the classical algorithm and the preconditioned algorithm, At = 
0.01, Ax = 1 /2048, Ay = 1 /64_ 


N 


256 

512 

Classical algorithm 

m 

Tnopc 

8 

417.9 

12 

344.9 

Preconditioned algorithm 

m 

Tpc 

5 

259.9 

6 

268.7 


We finish this subsection by some numerical tests for a time dependent potential V = 5(x 2 + y 2 )( 1 + 
cos(47rt)). We get similar conclusion by the results shown in table [9j 

In conclusion, the preconditioner allows to reduce significantly the number of iterations and the com¬ 
putation time. 
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Table 9: Number of iterations of the first time step and total computation time in seconds of the classical 
algorithm and the preconditioned algorithm, T = 0.5, At = 0.01, Ax = 1/128, Ay = 1/8. 


N 

2 

4 

8 

16 

32 

-^nopo P = 10 

31 

31 

31 

31 

31 

Ypc, p = -10 

9 

9 

9 

9 

9 

rp ref 

263.2 

^riopc 

272.6 

138.5 

72.4 

40.5 

19.5 

Tpc 

205.6 

101.4 

52.2 

29.4 

14.8 


7 Conclusion 

We applied the SWR method and the DDS method to the two dimensional linear Schrodinger equation 
with general potential. We proposed a new algorithm if the potential is a constant and a preconditioned 
algorithm for a general linear potential, which allows to reduce the number of iterations and the computa¬ 
tion time compared with the classical one. According to the numerical tests, the preconditioned algorithm 
is not sensitive to the transmission conditions (Robin, 5'™ de ) and the parameters in these conditions. 
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